PROGRAM MAPS
USE POLYMORPHIC_COMPLEXTAYLOR
TYPE(DAMAP) TOTAL_DA_MAP,partially_inverted_map
REAL(DP) ANGLE
type(pbfield) h
type(taylor) g
integer, allocatable :: j(:)
integer mf,me
mf=20

open(unit=mf,file='results.txt')

CALL INIT(NO1=3,ND1=1,NP1=0,NDPT1 =0)     !   <------------------ init for maps in ND1 degrees of freedom

CALL ALLOC(TOTAL_DA_MAP, partially_inverted_map)

ANGLE=acos(0.6d0)

call alloc(h); call alloc(g);

! This illustrates partial inversion by computing a generating function manually so to speak

! TOTAL_DA_MAP = Rot(30) exp(:x**3:) Identity 

TOTAL_DA_MAP=1
h= 1.d0.mono.'3'
TOTAL_DA_MAP=texp(h,TOTAL_DA_MAP)
h=(-angle/2.d0)*((1.d0.mono.'2')+(1.d0.mono.'02')) 
TOTAL_DA_MAP=texp(h,TOTAL_DA_MAP)



allocate(j(c_%nd2))
j=0
j(2)=1

 partially_inverted_map=TOTAL_DA_MAP**j

write(mf,*) "  "
write(mf,*) " The map ==> Rot(30)exp(:x**3:)Identity "
call print(TOTAL_DA_MAP,mf) 

write(mf,*) "  "
write(mf,*) " Partially inverted map "
call print(partially_inverted_map,mf) 


write(mf,*) "  "
write(mf,*) " The generating function G(q,p^final) "
! The subroutine intd_taylor extracts a generating function in (q,p^final) or (q^final,p) if factor=1.d0
! If factor = -1.d0, it extracts a Lie Polynomial
! The input is an array of Taylors (Not a DAMAP !!! )
call intd_taylor(partially_inverted_map%v,g,factor=1.d0)

call print(g,mf) 

deallocate(j)
call KILL(h);
call KILL(g);
CALL KILL(TOTAL_DA_MAP, partially_inverted_map)
close(mf)
END PROGRAM MAPS

